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The Bak-Sneppen model is a simple stochastic model of evolution that exhibits self-organized criticality and 
for which few analytical results have been established. In the original Bak-Sneppen model and many subse- 
quent variants, interactions among the evolving species are tied to a specified topology. We report a surprising 
connection between Bak-Sneppen type models and more tractable Markov processes that evolve without any 
reference to an underlying topology. Specifically, we show that in the case of a large number of species, the 
long time behaviour of the fitness profile in the anisotropic Bak-Sneppen model can be replicated by a model 
with a purely rank-based update rule whose asymptotics can be studied rigorously. 
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In (T), Bak and Sneppen introduced a very fruitful and 
simple model of evolution that exhibits interesting dy- 
namics but has proved surprisingly hard to analyse. The 
classical Bak-Sneppen (BS) model is a simple stochastic 
coarse-grained model of evolution of an ecosystem con- 
sisting of a fixed number N of evolutionary niches organ- 
ised in a ring. Each niche is occupied by a species with 
a particular^fness value in [0, 1]. Direct inter-species in- 
teractions (predation, competition, etc.) occur only be- 
tween species in neighbouring niches. The dynamics 
of the system is driven by the removal (extinction) of 
the least fit species in the entire system, whose niche is 
taken over by a new species; the extinction of the least fit 
species induces changes in the fitnesses of the species in 
the two neighbouring niches. In this letter, we show that 
a process whose update rule is defined solely in terms of 
the ranks of the fitness values, without any reference to 
topology of interactions, exhibits asymptotic behaviour 
similar to the (anisotropic version of the) BS model as 
well as self-organized criticality [2|. We call processes of 
this type rank-driven processes (RDPs) and analyse them 
rigorously in [ 3 1 . RDPs are of independent mathemati- 
cal interest and can be used to define new evolution mod- 
els. Despite the considerable impact of the Bak-Sneppen 
model on the physics community and beyond H, so far 
only a small number of rigorous results on the BS model 
have been obtained, such as those of Meester and Zna- 
menski [5 1 on the non-triviality of the steady-state distri- 
bution. In this letter, by exploiting the tools for analysis 
of RDP models developed in 0, we provide an approach 
for establishing new results in this active area. 

The BS model [1] is a discrete-time Markov process 
which advances every time there is a species extinction 
event. Each species occupying the N niches is initially 



assigned a fitness Xk G [0, 1] (fc G {1, . . . , N}) cho- 
sen independently from the uniform distribution on the 
unit interval, U[0, 1]. At each step of the algorithm, we 
choose the smallest of all the Xf., Xkmin sav > an d replace 
Xkmin and its two nearest neighbours Xkmin±i (indices 
calculated modulo N) by new independent U[0, 1] ran- 
dom numbers. In simulations with large N, the marginal 
distribution of the fitness at any particular niche is seen 
to evolve to a U[s* , 1] distribution, with s* s» 0.667. 

A number of variants of Bak and Sneppen's original 
model have been introduced which evolve according to 
different topological criteria. One simple variant is the 
anisotropic Bak-Sneppen (aBS) model, in which, in ad- 
dition to the least fit species, only its right-hand nearest 
neighbour is replaced. The aBS model is the main focus 
of this letter because while it simplifies calculations, it 
preserves the key qualitative phenomena of the original 
BS model. For example, the aBS model also gives rise 
(according to large-A^ simulations) to a threshold value 
s* sa 0.724 (6]. One contribution of the present letter 
is to propose a characterization for s* in terms of osten- 
sibly simpler quantities associated with aBS. Our argu- 
ments connecting rank-driven processes to aBS apply to 
BS too. 

Another variant on the BS model which eliminates 
topology is the mean-field version analysed in CH91, in 
which one replaces the smallest fitness and K — 1 ran- 
domly chosen other ones. Below we show that such mod- 
els fall within the RDP framework. 

Consider a process in which at each update the species 
with the smallest fitness and the i?-th ranked fitness are 
replaced, where R is a random variable on {2, 3, ... , N} 
(R = N corresponding to the largest fitness) sampled 
independently from a distribution P[R = k] = /jv(&) 



where /jv(fc) > and Y2k=2 flf(k) = 1- This is 
an example of a rank-driven process. The complexity 
of this RDP is intermediate between that of the aBS 
model and the mean-field model of fl7r|9); the latter is 
the special case of a RDP with /j\r(fc) = jrzj for all 
k G {2, . . . , TV}. The RDP has the advantage over aBS 
that it can be analysed rigorously. We have strong nu- 
merical evidence that for a judicious choice of /tv(/c) this 
simpler model can replicate the asymptotic behaviour of 
aBS. 

Specifically, one can choose /jv(^) to be fff s (k), the 
empirical distribution of the rank of the second chosen 
site in aBS: if we let P(k, M) be the number of times the 
k-th ranked element, k > 2, is the right neighbour of the 
smallest element in M iterations of the aBS algorithm, 



quantity is 



ff? s (k) = Jim jjP(k,M). 



(1) 



Heuristically, we expect that given suitable ergodicity 
properties for the aBS Markov process on the uncount- 
able state space [0,1]^, this limit will exist with proba- 
bility one. 

Little is known analytically about /j^ BS , due to the dif- 
ficulty of the aBS model, but /^ BS can be accurately 
numerically computed. Figure [TJ shows simulation es- 
timates of ffp s {k) for small values of k and different 
values of TV. 
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The RDP with = /^ BS computed numerically 
replicates many aspects of the dynamics of aBS. Specif- 
ically, rigorous results show that such an RDP exhibits 
a (large-TV, long-time) threshold at some s* which is 
an explicit function of f^: for the particular choice 
In = //v BS > using our numerical estimates for /j^ BS 
gives rise to a value for s* which is very close to the 
value computed directly from simulations of aBS. 

In O, by considering the random walk associated with 
an RDP defined in terms of f^, we show that a crucial 



a = lim lim 

n—>oo N— >■» 



(2) 



k=2 



assuming that the TV-limit exists. Here a <E [0, 1] mea- 
sures the "atomicity" of /jy as TV — > oo. For example, 
for the mean field aBS of [7-9 1 one has a = 0, while if 
we always replace the smallest and the second smallest 
elements, a = 1. The main result of [3 1 is that the thresh- 
old in the limiting (TV — > oo) stationary distribution of x 
values in the RDP is given by 

✓ -!±=. 

A second result of [ 3 1 shows that the limiting marginal 
distribution at stationarity is U[s* , 1], where s* is given 
by ([3]), provided that the selection distribution /jv is 
"eventually uniform" in the sense that 

Mfe) « ^ (4) 

for k sufficiently large. This condition is satisfied with 
a = for the mean-field aBS model, showing that the 
limiting distribution is indeed U[l/2, 1] in that case, as 
indicated by [8 1. 

From the results of Figure [T] we see that for a given 
TV, fff s {k) decays rapidly for small k before settling 
down to a uniform value. In fact, it appears that there is a 
constant C such that f^P s (k) = C/TV for large enough 
k. Thus the numerical evidence supports the eventual 
uniformity condition (|4]). Hence a = 1 — C. Numer- 
ical results give a ~ 0.445 and hence s* w 0.723, in 
close agreement with the simulations of |6|. Note that 
/^ BS (2) fa 0.209 for all the values of TV in Figure[l] 

Following [ 1 1 we define the length of an s-avalanche 
to be I if the number of consecutive steps for which the 
smallest fitness value stays below s is I. We compute the 
distribution n(£) of s-avalanche lengths for aBS and our 
RDP. Representative distributions are given in Figure [2] 
As s approaches s* we find that n(l) shows the power 
law behaviour characteristic of self-organized criticality, 
although there is a small but clear difference in the expo- 
nents of the two processes. 

The class of RDPs defined by /jy described above is 
contained in a wider class of processes [3|, which we 
now define. In this generality, an RDP is a discrete-time 
Markov process on the TV-simplex 

Aw = {0(i), ■ • -,x (N) ): < X( X ) < < x (N) < 1}; 

■ ■ ■ , %(N) are trie (increasing) order statistics of 
x\,..., x 7v ■ The RDP evolves according to the follow- 
ing Markovian rule. At each step, K of the Xk -values 
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FIG. 2. Size distribution n(£) of s avalanches in aBS (solid 
line) and RDP (dashed) for s = 0.5, 0.68 and N = 1000. 

are selected, according to rank, by sampling (without re- 
placement, according to some specified probability dis- 
tribution) from {1, 2, . . . , N}; that is, the sample from 
{1,2, ...,N} specifies the X(k) th at are chosen. The 
chosen Xfc-values are replaced by new (independent) 
U[0, 1] values. These processes are discussed in more 
detail in (3] . 

The numerical evidence reported above leads to sev- 
eral interesting problems for further investigation, not 
least the strong suggestion that the RDP with /jv = ff^ s 
given by ([T]i is closely related to aBS itself. The exact re- 
lationship of the two processes remains to be character- 
ized rigorously. If one wished to define a Markov process 
on Ajv whose stationary distribution coincided with the 
projection onto A^v of the stationary distribution of aBS, 
a natural candidate would be a RDP with state-dependent 
selection distribution: instead of a single /jv( - ) one 
would have a family /jv ( • ; x) of selection distributions 
conditioned on the state x G A at- Thus, assuming it ex- 
ists, one would take /jy( ■ ; x) to be fff* ( • ; x), the sta- 
tionary distribution for aBS of the right-neighbour of the 
smallest element conditional on the projection of the cur- 
rent state onto An being x. The fact that the numerical 
evidence described above suggests that one can proceed 
not with a state-dependent RDP based on ffp s ( • ; x) but 
with the simpler RDP based on ffp s ( • ) (which is an 
average of the f^p s ( ■ ; x)) seems to point to some im- 
portant underlying property of aBS itself. Two possible 
explanations are: 

(a) fjP S ( - ) = fpp S {- ; x) for all x, i.e., at station- 
arity there is some independence between the or- 
der statistics and the permutation that maps sites 
to ranks; or 

(b) f^p s {-',x) satisfies (uniformly in x) the same 
asymptotic conditions as /^r BS ( • ) that are central 
to the limit behaviour, namely Q and Q. 



The stronger fact (a) would suggest that the stationary 
distribution of the RDP coincides with the projection of 
the stationary distribution of aBS onto A at, so that the 
two processes share the same detailed equilibrium prop- 
erties. The weaker fact (b) would suffice to explain why 
the two processes share the same threshold and charac- 
teristic f/[s*,l] limit distribution. We remark that the 
distributions ffp s ( ■ ; x) seem to be very difficult to eval- 
uate numerically. 

Finally, we give a very brief indication of the origin of 
the threshold formula fij; see (3) for details. Consider 
the s-counting process N t (s) defined to be the number 
of Xfc-values in the interval [0, s] after t iterations of the 
RDP defined by /jy. Then N t (s) is a Markov chain on 
the finite state-space {0, 1, ... , N}. The threshold s* re- 
lates to the limiting (t — > oo then N — > oo) marginal 
distribution of an arbitrary Xk. The probability that a 
randomly chosen x^-value is less than s is E[N t (s)]/N 
(where E denotes expected value). Thus a natural way 
to define a threshold s* is 

s* = sup{s > : lim lira N- x E[NAb)\ = 0}; 

N—foo t->oo 

the t-limit exists by Markov chain limit theory and it can 
be shown that the TV-limit exists too, so that s* is well 
defined [3 1. 

To evaluate s*, we compute the mean drift of N t (s): 

E[N t+1 (s) - N t (s) | N t (s) = n] 

= 2s-(l + F N (n))l{n>0}, (5) 

where Fpj(n) = Ylk=2 /jv(fc)> 1 denotes an indicator 
function, and an empty sum is 0. Heuristically, for large 
N and large n, Fn(ii) ps a by (j2j) so that this drift is ap- 
proximately 2s — 1 — a, and setting this equal to zero 
gives (|3j. One expects that the drift being zero indi- 
cates the threshold behaviour, because a positive (neg- 
ative) drift would mean N t (s) increases (decreases). In 
this argument there are several limits involved (n, N, t 
all going to oo) that need to be handled with care. We 
exploit techniques from Markov process theory, such as 
Foster-Lyapunov conditions ifTUl . to do this: we refer to 
|3| for the details. 

In conclusion, we have indicated how the distribution 
/^r BS (fc) and the quantity a of (|2j capture the build-up of 
correlations in Bak-Sneppen type algorithms, the thresh- 
old behaviour of which can be analysed exactly by con- 
sidering the N t (s) Markov process on a countable state 
space. 

The class of RDPs that we have introduced is of inter- 
est in its own right. Numerical evidence suggests that by 
choosing as parameter for the RDP an appropriate statis- 
tic (/^P ( • )) of aBS, one can replicate the asymptotic 
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behaviour of aBS by the RDP, for which one can prove 
rigorous results more easily. The remaining analytical 
challenge is to clarify the relationship between aBS and 
the RDP. This involves at least two main parts: (i) prov- 
ing the existence of the distributions /JJ, BS given by (jlj 
and of the limit a defined by |2]); and (ii) determining the 
property of aBS that allows us to use fff s ( ■ ) instead 
of the conditional version /JJr BS ( • ;x). If one can make 
precise the connection between aBS and the RDP, one 
should be able to transfer rigorous results for RDPs [ 3 1 
to aBS. In respect to challenge (i) above, it is interest- 
ing to note that if an explicit description of ffp s could 
be obtained, one might be able to obtain an explicit for- 
mula for the threshold s* via Q and ([3]). Finally, we note 
that in the case of the classical BS process the same ar- 
guments apply, though now /at is a function of two vari- 
ables, f N (k,£), k € {2, ... N-l}Je {k+1, ...,N}. 

MG would like to acknowledge fruitful discussions 
with Gregory Berkolaiko, Jack Carr, Oliver Penrose, and 
Michael Wilkinson. 
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